Parametrization of the Hybrid Potential for Pairs of Neutral Atoms 
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The hybrid form is a combination of the Rydberg potential and the London inverse-sixth-power 
energy. It is accurate at all relevant distance scales and simple enough for use in all-atom simulations 
of biomolecules. One may compute the parameters of the hybrid potential for the ground state of a 
pair of neutral atoms from their internuclear separation, the depth and curvature of their potential 
at its minimum, and from their van der Waals coefficient of dispersion C% . 

PACS numbers: 34.20.Cf (Interatomic potentials and forces), 87.15.-v, 82.35.Pq, 77.84.Jd 



I. THE HYBRID FORM 

Commercial molecular-modeling codes use pair poten- 
tials chosen more for speed than for accuracy. They use 
the harmonic potential 



V h (r) 



-Eo + - (r - r ) 2 



(1) 



for covalcntly bonded pairs of neutral atoms and the 
Lennard-Jones potential [l[ 



V LJ (r) = E 



(2) 



for unbonded pairs. The recently introduced hybrid po- 
tential @ 



V(r) 



- b r 



(i 



Ca 



dr 



(3) 



is nearly as fast and much more accurate. It is fast 
enough for use in all-atom simulations of biomolecules 
and accurate at all biologically relevant distance scales 
when its parameters are suitably chosen. How to choose 
them is the focus of this paper. Formulas are derived 
for a, b, c, & d in terms of the internuclear separation 
ro, the depth Eq and curvature k of the potential at 
its minimum, and the van der Waals coefficient Cq of 
the pair. The hybrid potential therefore is applicable to 
pairs of neutral atoms for which no empirical potential 
is available. Its eventual incorporation into TINKER Q, 
Amber Q , and other such codes is a goal of this work. 

The hybrid potential is a combination of the Rydberg 
formula used in spectroscopy and the London formula 
for pairs of atoms. The terms involving a, b, and c were 
proposed by Rydberg to incorporate spectroscopic data, 
but were largely ignored until recently. The constant Cq 
is the coefficient of the London tail. The new term dr~ 6 
cures the London singularity. As r — > 0, V(r) — > a, 
finite; as r — > oo, V(r) approaches the London term, 
V(r) —* — Ce/r 6 . In a perturbative analysis @, the 



a, b, c terms arise in first order, and the Cq term in second 
order. 

When fitted to spectroscopically determined potentials 
for the ground states of H2, N 2 , O2, Ar-Ar, and Kr-Kr, 
the hybrid form is four orders of magnitude more ac- 
curate than the harmonic and Lennard-Jones potentials 
and five times more accurate than the Morse J6( , Varn- 
shi Q , and Hulburt-Hirschfelder [j| potentials 0] . It also 
yields accurate second virial coefficients and heats of va- 
porization Q . Its simplicity recommends it as a teaching 
tool and as a practical form for computation. 

How does one find the parameters a, b, c, & d when an 
empirical potential is not available? For many pairs of 
neutral atoms, spectroscopists have measured the inter- 
nuclear separation ro, the well depth Eq = |V(ro)| and 
curvature k = V (ro) at the minimum of the potential, 
and the London coefficient Cq. These input parameters 
are discussed in Sec. |TT] with an emphasis on the curva- 
ture k and its relation to the vibrational frequency of 
the ground state and to the energy Dq needed to dis- 
sociate the ground state. Section IIIII develops a rule of 
thumb for the parameter d to which the hybrid form is 
relatively insensitive. Section ITVl derives formulas for the 
hybrid parameters a, b, & c in terms of ro, Eq, k, Cq, & d 
and plots the resulting hybrid potentials for the 11 pairs 
H 2 , N 2 , 2 , NO, OH, I a , Li 2 , Na 2 , K 2 , Ar-Ar, & Kr-Kr. 



II. INPUT PARAMETERS 

In addition to the well depth Eq = \V(ro)\, the inter- 
nuclear separation ro at the minimum, and the London 
or van der Waals coefficient Cq, the curvature k = V (r ) 
of the potential at ro often is available either directly or 
in the guise of the energy of dissociation Dq . 

Near ro , the motion of the internuclear separation r is 
described by the lagrangian 



■ /ir' 



-^fc(r-ro) 2 



(4) 



in which fi = mim 2 /(mi + m 2 ) is the reduced mass of 
the two atoms of mass mi and m 2 . To lowest order in 
r — ro, the ground state then has energy 
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in which u> is the (angular) vibrational frequency of the 
ground state. The curvature fc is related to u> by 



2 

fiU) . 



(6) 



In the chemical-physics literature, the depth En and 
separation ro arc labeled D e and r e , and the angular 
frequency u> = 2ttv is expressed as a frequency u e in 
inverse centimeters: u> e — lu/(2ttc) where c is the speed 
of light. In these terms and in cgs units, k = Ait 2 c 2 ji us 2 . 
The energy of dissociation is Do = — E g = En — fkv/2, and 
so the curvature k — V (ro) is related to the difference 
between it and the well depth En = D e by 



k = 4 M (£ - Dnf/h 2 = A f i(D e - Dnf/h 2 



(7) 



The values of the input parameters are listed in Ta- 
ble |TJ The internuclear separation r at the minimum 
of the potential well V(r) and its depth E = V(r ) are 
from [| for H 2 , NO, OH, & I 2 , and from [3 for 2 . 
The curvature k = V (ro) of the potential at its min- 
imum is from for H 2 , 2 , NO, OH, & I 2 . The van 
der Waals coefficient Cg is from fill for H 2 ; from [l2| for 
NO; from [3 for OH; and from [l| for 2 , I 2 , Ar-Ar, & 
Kr-Kr. The values of rg. En, & k for Ar-Ar and Kr-Kr 
respectively are from [ll| and The values of rn En, 
fc, fc Cr for N 2 , Li 2j Na 2 , & K? respectively are from fl7j ]; 
[H; [11, [H], & [U; and gg, |U, & Q. 



III. A RULE OF THUMB 

The d-term of the hybrid form is a trick to avoid the 
1/r 6 singularity of the London term. It has no other 
justification. Luckily, the quality of the fit is not sensitive 
to the precise value of the parameter d, and so we need 
only a good rule of thumb for that parameter. To develop 
an approximate formula that estimates this parameter, 
we first find values of a, b, c, & d that provide good fits 
of the hybrid form V(r) to empirical potentials for the 
11 pairs H 2 , N 2 , 2 , NO, OH, I 2 , Li 2 , Na 2 , K 2 , Ar-Ar, 
& Kr-Kr of neutral atoms. 

Empirical potentials obtained from spectroscopic 
data [HI [H Q by the RKR (Rydberg [H, Klein [26[, 
method are available from 0] for H 2 , NO, & 



Rees [27 
I 2 ; from 
from [15 
from [15 
Section 



281 ] for OH; from 
for Li 2 ; from [19 
for Ar-Ar; and from [l6 
IVI contains formulas 



17J for N 2 ; from llOU 
for Na 2 ; from [29| 
for Kr-Kr. 

cna rrn m 



for 2 ; 
for K 2 ; 

& HO]) 

for a, b, &c in terms of En, rn, k, Cg, &c? that ensure 
that the hybrid form V(r) goes through the minimum 
(ro,— En) with curvature fc. For each pair of atoms, I 
found the value of d that best fits the hybrid form V(r) 
to the empirical potential of the pair for r > f rn in 
which < / < 1. The resulting values of a, b, c, & d 
are listed in Table [Til along with the fractions / and the 
root-mean-square (rms) errors AV(r) for r > f rn. Fig- 
ures [THTTI show that the fitted hybrid forms Vt (dashes, 
dark blue) nicely follow the points (diamonds, cyan) of 



TABLE I: The depth E = V(r ) of the potential well V(r), 
the internuclear separation ro and the curvature k = V"(ro) 
at the minimum r = ro, and the van der Waals coefficient 
for 11 pairs of neutral atoms. 







ro (A) 




t-6 lev A ) 


TT 




7/117 


OO.ooDl 


O OQOQO 

o.ooooo 


14 N 2 


9.8995 


1.09768" 


143.2245 


14.382 


o 2 


5.2136 


1.2075 


73.4726 


9.3215 


NO 


6.609 


1.1590 


99.5676 


11.245 


OH 


4.624 


0.9707 


48.6196 


6.854 


h 


1.5571 


2.668 


10.7378 


230.05 


'Li 2 


1.0559 


2.6730 


1.5752 


829.33 


Na 2 


0.74664 


3.0786 


1.0706 


929.76 


K 2 


0.55183 


3.9243 


0.61375 


2328.6 


Ar 2 


0.01234 


3.757 


0.0691 


38.4213 


Kr 2 


0.01735 


4.017 


0.0896 


77.6791 



"Extra digits are included to avoid round-off errors. 



TABLE II: The values of the coefficients a, b, c, & d that fit 
the hybrid potential V{r) (Eq. (J3j) ) to RKR data for H 2 , N 2 , 
2 , NO, OH, l 2 , Li 2 , Na 2 , & K 2 and to empirical potentials 
for Ar-Ar & Kr-Kr while respecting Eqs. ([Tol [171 IT91 fe l2u) . 
The rms errors AV(r) for r > fro are also listed. 





a (eV) 


b (A-) 


c (A- 1 ) 


d (A 12 ) 


f, AV(eV) 


H 2 


47.796 


2.9632 


2.5406 


12.2 


0.68, 0.086 


N 2 


3752.644 


4.3533 


1.1777 


34.8 


0.81, 0.088 


o 2 


2901.580 


4.2173 


1.0510 


59.8 


0.84, 0.039 


NO 


3809.497 


4.4196 


1.0943 


47.0 


0.81, 0.077 


OH 


377.804 


3.6909 


1.4668 


32.5 


0.80, 0.023 


l 2 


14361.15 


2.8013 


0.4351 


2.08e5" 


0.88, 0.039 


Li 2 


199.481 


1.6200 


0.5101 


2.85e6 


0.86, 0.016 


Na 2 


231.900 


1.5311 


0.4292 


9.40e6 


0.81, 0.009 


K 2 


325.051 


1.3409 


0.3269 


9.94e7 


0.85, 0.008 


Ar 2 


1987.943 


2.6517 


0.2978 


8.10e7 


0.70, 0.00009 


Kr 2 


1875.462 


2.3661 


0.2789 


5.72e8 


0.90, 0.00018 



6 2.083e5 = 2.083x10°. 



the empirical potentials for the 11 pairs of neutral atoms 
H 2 , N 2 , 2 , NO, OH, I 2 , Li 2 , Na 2 , K 2 , Ar-Ar, & Kr- 
Kr. The fitted Vf(r)'s go through the empirical minima 
with the right curvatures and closely trace the empirical 
potentials, at least for r > frn. Figure [1] for molecu- 
lar hydrogen adds the harmonic potential V/,. (dot-dash, 
magenta), which fits only near the minimum at r = r . 
Figure [2] for molecular nitrogen and Fig. [TD] for a pair 
of argon atoms include the Lennard-Jones potential Vlj 
(dot-dash, magenta), which is accurate only near ro- 

The fitted values of the parameter d in Table [TTJ arc 
less well defined than those of a, b, and c. The empirical 
rule of thumb 



7.1 A 1 



2.89- 6 



0.468-^-^ 8 
ev A 



gives the 11 values of d listed in Table ITjTI and roughly 
approximates those of Tabic [TTJ 
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FIG. 1: The hybrid form V with the calculated coefficients (a, 
b, c, & d, Table IIII|I (solid, red) fits the RKR spectral points 
for the ground state of molecular hydrogen (diamonds, cyan) 
and gives the correct London tail for r > 3 A nearly as well as 
does the hybrid form Vf with the fitted coefficients (a, b, c, & 
d, Table HH (dashes, dark blue). In all the figures, the values 
of C*6 used in V and Vf are from Table U The harmonic form 
Vh (dot-dash, magenta) fits only near the minimum. 
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FIG. 2: The hybrid form V with the calculated coefficients 
(a, b, c, & d, Table IIII|I (solid, red) fits the RKR spectral 
points for the ground state of molecular nitrogen (diamonds, 
cyan) nearly as well as does the hybrid form Vf with the fitted 
coefficients (a, b, c, & d, Table W\\ (dashes, dark blue). The 
Lennard- Jones form Vlj (dot-dash, magenta) fits only near 
the minimum. 



TABLE III: Values of the coefficients a, 6, c, and d obtained 
by Eqs. (JSJ US] ED EH horn the values of E , r , k, and 

C 6 of Table U 





d i^ev ) 


h (k- 1 \ 
b (A ) 


C (A J 


a (A ) 


f AT/ 
I, ZA V 


TT 

ri2 


40. Ul 


o ori7 
z.yu / 


z.ODOo 


ID. / 


U.Do, U.Uo/ 


IN 2 


,int;Q no 




1 1 7ftO 
1.1 / DZ 




n si n 1 1 1 


o 2 


2868.48 


4.246 


1.0539 


40.7 


0.84, 0.068 


NO 


4040.42 


4.496 


1.0946 


34.0 


0.81, 0.156 


OH 


491.96 


3.942 


1.4478 


19.8 


0.81, 0.053 


h 


16125.67 


2.832 


0.4350 


2.79e5 


0.88, 0.041 




148.93 


1.516 


0.5161 


4.01e6 


0.86, 0.028 


Na 2 


278.20 


1.595 


0.4275 


7.33e6 


0.81, 0.012 


K 2 


359.65 


1.381 


0.3273 


7.30e7 


0.85, 0.012 


Ar 2 


4994.79 


2.921 


0.2959 


3.12e7 


0.70, 0.012 


Kr 2 


9610.07 


2.805 


0.2759 


6.23e7 


0.90, 0.0004 
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FIG. 3: The hybrid form V with the calculated coefficients 
(a, 6, c, & d, Table IIII[l (solid, red) fits the RKR spectral 
points for the ground state of molecular oxygen (diamonds, 
cyan) nearly as well as does the hybrid form Vf with the fitted 
coefficients (a, b, c, & d, Table HT) (dashes, dark blue). 



IV. THREE FORMULAS 

Suppose we use the rule (0) for the fudge-factor d and 
take the London coefficient C§ from Table U How do 
we find the parameters a, b, and c that make the hybrid 
potential V(r) have its minimum at r = ?*o with V(ro) = 
— Eq and with curvature V (to) = fc? 

Let us write the hybrid form as the sum 

V(r) = v(r) + w(r) (9) 
of the a, b, c terms v(r) 

v(r) = ae- br {l - cr) (10) 
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FIG. 4: The hybrid form V with the calculated coefficients (a, 
b, c, & d, Table IIII|I (solid, red) fits the RKR spectral points 
for the ground state of nitric oxide (diamonds, cyan) nearly 
as well as does the hybrid form Vf with the fitted coefficients 
(a, b, c, & d, Table HI)) (dashes, dark blue). 
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FIG. 6: The hybrid form V with the calculated coefficients 
(a, 6, c, & d, Table IIII[) (solid, red) fits the RKR spectral 
points for the ground state of molecular iodine (diamonds, 
cyan) nearly as well as does the hybrid form Vf with the 
fitted coefficients (a, b, c, & d, Table ITT)) (dashes, dark blue). 
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FIG. 5: The hybrid form V with the calculated coefficients 
(a, 6, c, & d, Table IIII)l (solid, red) fits the RKR spectral 
points for the ground state of the hydroxyl radical (diamonds, 
cyan) slightly less well than the hybrid form Vf with the fitted 
coefficients (a, b, c, & d, Table HT)| (dashes, dark blue). 




FIG. 7: The hybrid form V with the calculated coefficients 
(a, b, c, & d, Table IIII|) (solid, red) fits the RKR spectral 
points for the ground state of the lithium dimer (diamonds, 
cyan) nearly as well as does the hybrid form Vf with the fitted 
coefficients (a, b, c, & d, Table Hi)) (dashes, dark blue). 



and the C§, d terms w(r) 



dr~ 6 ' 



(11) 



Since Cq and d are given, the function w(r) and its deriva- 



tives 



w '(r) = 6C 6 r 5 



„12 



d 



(r 12 + df 



(12) 
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FIG. 8: The hybrid form V with the calculated coefficients (a, 
b, c, & d, Table IIII|I (solid, red) fits the RKR spectral points 
for the ground state of the sodium dimer (diamonds, cyan) 
well, but not quite as well as does the hybrid form Vf with 
the fitted coefficients (a, b, c, & d, Table [TTJ) (dashes, dark 
blue). 
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FIG. 9: The hybrid form V with the calculated coefficients (a, 
b, c, & d, Table ITTT|) (solid, red) fits the RKR spectral points 
for the ground state of the potassium dimer (diamonds, cyan) 
well, but not quite as well as does the hybrid form Vf with 
the fitted coefficients (a, b, c, & d, Table [TTJ) (dashes, dark 
blue). 
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and 



w"(r) 



-6C 6 r 



4 (7r 12 ~d){r 12 -5d) 



(r 12 + d) 3 

are determined. The condition that V(rg) = - 
is 

ae- br "(l-cr Q ) = -w(r )-E 
which implies that a is 

a = -(w(r Q )+E a )e br °/(l~cr ). 
The condition that V'(ro) = is 

ae- br °{b+c-bcr )=w'{r Q ) 
which together with (fT5|) gives b as 

b = -c/(l - cr ) - «/(r )/(«>(ro) + 
Finally, the condition V (fo) = fc is 

a be- bro (b-bcr +2c) = k-w" (r ) 



(13) 
-Eo then 

(14) 
(15) 
(16) 
(17) 
(18) 



which with (|15|) for a and p7p for b is a quadratic equa- 
tion for c with roots 



(w(r ) + E o y 



w'(r ) 2 + (k- w"(r ))(w(r ) + E Q ) 
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FIG. 10: The hybrid form V with the calculated coefficients 
(a, 6, c, & d, Table IHIfl (solid, red) fits the RKR spectral 
points for the ground state of the argon dimer (diamonds, 
cyan) nearly as well as does the hybrid form Vf with the fitted 
coefficients (a, b, c, & d, Table HT|) (dashes, dark blue). The 
Lennard- Jones form Vlj fitted to the minimum (dot-dash, 
magenta) is too low for r > 4 A. 
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FIG. 11: The hybrid form V with the calculated coefficients 
(a, b, c, & d, Table (solid, red) fits the RKR spectral 
points for the ground state of the krypton dimer (diamonds, 
cyan) nearly as well as does the hybrid form V/ with the fitted 
coefficients (a, b, c, & d, Table HT|| (dashes, dark blue). 



the depth and curvature of their potential at its mini- 
mum, and from their van der Waals coefficient C§. The 
hybrid potential therefore is applicable to pairs of neu- 
tral atoms for which no empirical potential is available. 
Given the differences between it and the Lennard-Jones, 
harmonic, Morse, Varnshi, and Hulburt-Hirschfeldcr po- 
tentials, it would be worthwhile to examine the conse- 
quences of these differences in Monte Carlo searches for 
low-energy states of biomolecules and in numerical sim- 
ulations of phase transitions and reactions far from equi- 
librium. 



Acknowledgments 

Thanks to G. Groenenboom for sending me a copy of 
the potential of [n| and to S. Atlas, S. T. P. Boyd, H.-B. 
Broeker, D. Cromer, K. Dill, M. Fleharty, S. C. Foster, 
R. J. Le Roy, W. J. Meath, V. A. Parsegian, R. Pas- 
tor, R. Podgornik, E. Riley, B. Rivers, D. Sergatskov, 
J. Thomas, S. Valone and W. Zemke for advice. 

APPENDIX: THE CHOICE OF SIGN 



The minus sign implies 

c> l/r 1 -era < (20) 
while the plus sign implies 

c < l/r <S=> 1 - cr > 0. (21) 

A covalently bonded pair will have v(ro) < 0, which 
implies 1 — c?* < and so c > l/r , which entails the 
minus sign. One may choose the minus sign for the 11 
pairs of neutral atoms considered in this paper. A further 
discussion of the choice of sign appears in the Appendix. 

If one uses the rule of thumb (J3J) for d and the values 
of ro, Eq, k, & Ce from Table CO then one may find c from 
condition (TIT)]) with a choice of sign, and then b from (jTTJ) , 
and then a from (|15() . The resulting values of a, 6, c, & 
d for the 11 pairs of neutral atoms H 2 , N 2 , O2, NO, OH, 
I 2 , Li 2 , Na 2 , K 2 , Ar-Ar, & Kr-Kr are listed in Table ITnl 
for the case of a minus sign in Eq. p~9]) ; also shown are 
the rms errors AV(r) for r > fro for the f's of Table HT1 
Figures [TlfTTI show the hybrid form V with these values 
of a, b, c, & d and with the Cq's of Table |T] as solid 
red curves. They closely trace the empirical potentials 
(diamonds, cyan) and the fitted hybrid forms Vf (dashes, 
dark blue). 

V. CONCLUSIONS 

By using Eqs. (8j Q31 OH [HI & EH , one may build an 
accurate hybrid potential (J3j) for the ground state of a 
pair of neutral atoms from their internuclear separation, 



The choice of the minus sign in Eq. pU]) implies that 
the fudge factor d must exceed a certain lower limit. The 
parameter a is the value of the potential V(r) at r = 
and so must be positive. The minus-sign inequalities (|20j) 
imply that 1 — cro < 0. Thus Eq. (|15p will give a positive 
value for a only if 

w{r o ) + E o >0 (A.l) 
which by pip implies the lower limit 



on the fudge factor d. The values of d listed in Tables [TT1 
fc IIIII satisfy this constraint. 

The choice of the plus sign in Eq. (fl~5)) implies two upper 
limits on the fudge factor d. The plus-sign inequalities 
(|2Tj) imply that 1 - cr > 0. Thus Eq.([T5]) will give a 
positive value for a only if 

W (r ) + £ <0 (A.3) 

which by pip implies the upper limit 

d < - r 12 (A.4) 

on the fudge factor d. The right-hand-side of this in- 
equality is negative for the pairs N 2 , 2 , NO, & I 2 and 
less than unity for the pairs H 2 & OH, as may be seen 
from the values of Cq, r , & E listed in Table CO Thus 
the minus sign is required for the pairs N 2 , 2 , NO, & I 2 
and strongly indicated for the pairs H 2 & OH. 
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TABLE IV: The values of the coefficients a, b, c, & d that best 
fit the hybrid potential V(r) (Eq. ©) to RKR data for Li 2 , 
Na2, & K 2 and to empirical potentials for Ar-Ar & Kr-Kr 
while respecting Eqs. (|15l 1171 & I19p with a plus sign in (|19jl . 
AV(r) is the rms error in V for r > f ro. 





a (eV) 


b (A -1 ) 


c (A" 1 ) 


d (A 12 ) 


f, AV(eV) 


Li 2 


1136.21 


1.8218 


0.3225 


869. 


0.77, 0.037 


Na 2 


648.22 


1.5692 


0.3034 


2692. 


0.90, 0.004 


K 2 


500.31 


1.2673 


0.2486 


33770. 


0.70, 0.012 


Ar 2 


1754.01 


2.7054 


0.2620 


1.57e5 


0.70, 0.001 


Kr 2 


2725.43 


2.5622 


0.2485 


1.0e6 


0.90, 0.0003 



The second upper limit on d arises because both b and 



c must be positive (6 > because v(r) must vanish as 
r — > 00, and c > because c < implies u(r) > 0). Thus 
if the positive sign is chosen in the formula (|19|) for c, then 
1 — cr > 0, and also w(ro) + £0 < by (|A.3|) . and so 
Eq. (fT7| will give a positive value for b only if w'(r ) > 0, 
which by (|12|) implies the upper limit 

d<r* 2 . (A.5) 

The curves displayed in Figs. ITT-fTll all correspond to 
the choice of a minus sign in Eq. (Tli?|) . But for the five 
pairs L12, Na2, K2, Ar-Ar, & Kr-Kr, one also may get 
good fits to the empirical data by using the plus sign (|2l"j) . 
The resulting values of a, b, c, & d appear in Table IIVI 
The cfs obey ([Alii fe 1X5]). 
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